Low-energy dynamics of the two-dimensional 
S=l/2 Heisenberg antiferromagnet on percolating clusters 
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We investigate the quantum dynamics of site diluted S — 1/2 Heisenberg antiferromagnetic clus- 
ters at the percolation threshold. We use Lanczos diagonalization to calculate the lowest excitation 
gap A and, to reach larger sizes, study an upper bound for A obtained from sum rules involving 
the staggered structure factor and susceptibility, which we evaluate by quantum Monte Carlo sim- 
ulations. Scaling the gap distribution with the cluster length L, A ~ L~ z , we obtain a dynamic 
exponent z « 2D/, where Df = 91/48 is the fractal dimensionality of the percolating cluster. This 
is in contrast to previous expectations of z = Df. We argue that the low-energy excitations are due 
to weakly coupled effective moments formed due to local imbalance in sublattice occupation. 

PACS numbers: 75.40.Gb, 75.10.Jm, 75.10.Nr, 75.40.Mg 
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Antiferromagnets under doping with spinless vacancies 
exhibit interesting phenomena resulting from quantum 
fluctuations building on top of fluctuations of classical 
percolation. Various calculations for the two-dimensional 
(2D) Heisenberg model indicated that the quantum fluc- 
tuations for S = 1/2 arc strong enough to cause a 
quantum phase transition at some vacancy concentra- 
tiona = p c less than the classical percolation threshold 
p* pj. However, recent quantum Monte Carlo (QMC) 
simulations have demonstrated that long-range order sur- 
vives for all p < p* 0,0. Experiments on magnetically 
diluted layered cuprates agree with this, showing a low- 
temperature divergent correlation length up to the per- 
colation point 0. The fractal percolating cluster at p* 
is ordered Q. The dilution-driven phase transition is 
therefore a classical percolation transition of the static 
spin order, but quantum fluctuations present in the per- 
colating cluster lead to changes in quantities related to 
the low-energy dynamics. The critical exponents depend 
on classical percolation exponents and the dynamic ex- 
ponent z of the clusters p| , which thus becomes the focus 
of attention — it is the subject of this Letter. 

The singlet-triplet gap of a clean Z3-dimensional anti- 
ferromagnet on a lattice of even length L scales as L~ z 
with z — D The lowest excitations are the quan- 
tum rotor states responsible for symmetry breaking when 
L — > 00. Vojta and Schmalian recently generalized this 
result to percolating clusters; z = D f HI, where in 2D 
the fractal dimensionality Df — 91/48 0- QMC calcu- 
lations by Yu et al. || seem to confirm z = Df at p*, 
however in an indirect way from the finite-temperature 
scaling £ ~ T~ x l z of the correlation length. This form 
is expected for a quantum-critical system |9j but is not 
necessarily related to the finite-size gap of ordered clus- 
ters The exponent governing the finite-size gap is 
the one relevant to the T — dilution transition. Here 
we will take a more direct approach, studying the scaling 
of the excitation gaps of percolating clusters. Our main 
result of this study is that the dynamic exponent is much 
larger than expected; z w 2Df . We will argue that the 
lowest excitations are due to weak interactions between 



moments formed by local sublattice asymmetry. 
The diluted Heisenberg hamiltonian is 
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where are nearest neighbors on a 2D square lattice 
and et = 0, 1 with probability p and I— p. We focus on the 
percolation threshold, p = p* w 0.407 [7J and study the 
properties of clusters grown on an infinite lattice, keeping 
clusters of a target size n (see Ref. 3] ) . We first consider 
the gap A = E\ —Eq between the singlet ground state and 
triplet first excited state of clusters with an equal number 
of spins in the two (A, B) sublattices; n = tia + tib, 
iiA = Tifl. We also discuss ha ^ ub , for which the ground 
state is an S — \ua — ^s|/2 multiplet, and find a subtle 
difference in the nature of the low-energy excitations. 

Using Lanczos diagonalization we study clusters 
with n up to 20. In order to reach larger sizes we use 
a QMC method (stochastic series expansion, SSE H2). 
Due to statistical errors, it is not possible to evaluate very 
small gaps directly using QMC. We therefore consider an 
upper bound A* for the gap, 



A* = 2S'(7r, 7r)/x(7r, 7r) > A, 



(2) 



obtained from well known sum rules for the static struc- 
ture factor S(q) and susceptibility x(q); 



P duS(q,u;) = S(q), 
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where S(q, a;) is the dynamic structure factor. For a 
single mode A* = A. As in the undiluted system [||, 
we expect A* for Heisenberg clusters to scale as the true 
gap; A* ~ L~ z ~ n~ z / Df . At the very least, z obtained 
from A* will be a lower bound for the dynamic exponent. 

In the SSE simulations we reach the ground state using 
a (3 = J/T doubling procedure 0, checking all clusters 
individually for T — > saturation of xi^i 7r )- Because of 
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FIG. 1: (color online) (a) Singlet-triplet gap distribution for 
clusters of size n = 16, 18, 20. Representative clusters corre- 
sponding to the two peaks are shown, (b) Distribution of the 
upper bound A* for the gap for n — 16, 32, 64. 



the large z, extremely low temperatures are required and 
we have therefore only gone up to n — 192, for which 
some clusters already require j3 = 2 17 (S(tt, it) saturates 
at significantly lower f3 [2j). The computational effort 



(CPU time and memory) scales as n(3 

In Fig. n we show the probability distributions of the 
exact gap A of clusters with n = 16, 18, and 20 obtained 
on the basis of 4 x 10 4 samples for each size. We also 
show results for the upper bound A* for n — 16,32, and 
64 obtained from 6x 10 3 clusters each. A common feature 
for the small clusters is a double-peaked gap distribution. 
In the upper panel we show two illustrative clusters cor- 
responding to the two peaks. Identifying pairs of spins 
likely to form singlets (giving a rough singlet-product ap- 
proximation of the ground state) we find that the clusters 
of the lower peak always have "dangling spins" not be- 
longing to a local singlet (for any reasonable pairing; the 
locations of the dangling spins are of course not unique), 
whereas dangling spins are absent for the upper peak. 
The large-gap peak vanishes as n grows; clearly because 
it is unlikely for a large cluster not to have any dangling 
spins. Dangling spins arise from local imbalance in sub- 
lattice occupation. They should act as weakly interacting 
moments, with effective couplings decreasing with sepa- 
ration. One should hence expect small gaps for clusters 
with dangling spins at large separation. 

The role of isolated spins in the formation of long-range 
order has been pointed out by Bray-Ali and Moore |13| . 
Although some spins can be very weakly coupled (effec- 
tively) to the rest of the cluster, correlations between 
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FIG. 2: (color online) Finite-size scaling of the average of the 
gap A and the upper bound A*. The slope of the line fitted 
to (A*) is z/Df = 1.99. 



them can be stronger than within the backbone of the 
cluster. Arbitrarily weakly coupled moments formed by 
groups of spins can also correlate over long distances and 
hence even a "sloppy" fractal cluster can order at T = 
(a more thorough examination of long-range order was 
presented in [lj])- Here we will argue that excitations 
involving localized effective moments are lower in energy 
than the quantum rotor states considered in previous dis- 
cussions of the dynamic exponent Q . 

Fig. shows the size dependence of the disorder av- 
erages (A) and (A*). A linear fit to (A*) gives a slope 
Z /D f = 1.99(2), or z = 3.77(4). Although the small 
systems for which we have calculated the exact gaps are 
not yet in the scaling regime, the available (A) points 
reproduce well the shape of the (A*) curve in the range 
where they overlap. Comparing the A and A* distribu- 
tions for n = 16 in Fig. ^ we also see that they are very 
similar — A* is only shifted up slightly relative to A. 

The large z implies that the lowest excitations of these 
clusters are not quantum rotor states (which would give 
z = Dt w 1.89). To investigate the nature of the excita- 
tions, we have calculated the local susceptibility; 
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which can be used to define a local gap 0] Aj = 1/2% j, 
as in Eq. ■ The magnitude of the local gaps of three 
representative clusters are illustrated in Fig. [3] In the 
left cluster a few isolated sites with small gaps can be 
distinguished, whereas the small gaps are more spread 
out in the one to the right. In the central one there are 
isolated as well as more extended moments. Thus not 
only dangling spins but also other groups of spins can 
form effective moments. In all cases, the small-gap re- 
gions appear to be localized. Regions with large gaps 
can be seen where it is apparent that the connectivity fa- 
vors local singlet formation. Such singlets can contribute 
to isolating the effective moments from each other. 

Assuming that the low-energy physics of the clusters is 
governed by interactions between effective localized mo- 
ments, we can apply concepts of the well-known strong- 
disorder renormalization group (RG) moment decimation 
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FIG. 3: (color online) Local gaps Aj for n = 128 clusters, 
shown as a color coded mapping of — ln(Aj) to [0, 1] (for each 
cluster separately). The smallest gaps are A m i n /J = 0.027 
(left), 0.014 (middle), and 0.012 (right). The largest gaps are 
A max /J = 0.9 (left and right) and 0.22 (middle). 



procedure 0, 0, 0] to analyze the low-energy behav- 
ior. The effective interactions between the moments can 
be ferromagnetic or antiferromagnetic, and we therefore 
have to consider two moments combining into either a 
smaller or a larger moment in each step |lil (includ- 
ing the case of two equal moments forming a singlet and 
dropping out). On a cluster with ha = «b the decima- 
tion would ultimately leave two equal moments with an 
effective antiferromagnetic interaction. The lowest exci- 
tation is the singlet-triplet excitation of this pair. A key 
question is whether the typical size of the renormalizcd 
moments grows with the cluster size or remains finite. If 
the moments remain finite the low-energy excitations are 
localized in the sense of involving a finite number of spins. 
However, the distance spanned by a pair of moments can 
grow with the cluster size, leading to decreasing effective 
interactions and thus smaller gaps. 

If the excitations indeed are localized moment-pair ex- 
citations, the distribution of the first gap (i.e., A, which 
is here estimated as A* or the smallest of the local gaps, 
Amin) should be described by extremal value statistics 
[20l l2lj of the local gaps. We thus consider the distri- 
bution of Aj, taking into account a possible scaling with 
the system size, as suggested by the discussion above. 
We define a = AiL a ~ Ain a / D f. As shown in Fig. El for 
a rj 2.8 the distributions of e$ for different cluster sizes 
indeed collapse onto a single universal curve. For small 
ej, P(ei) ~ Aef , where uj w 1 and A is a constant. We 
want the distribution Pm (emin) of the smallest local gap, 
e mi „ = A min L a = min{ei,e 2 , • • • ,e M ], for large M. If 

Using 



A mi „L a = min{ei,e 2 , 
the excitations are localized we have M ~ L Df 
the probability of a scaled local gap e' < e, 

A 



P<(e) 



P(e')de' 



1 



(6) 



we find the typical e m j n from p<(e m in) ~ M . Then 
using M ~ L D f we get A min ~ L~ a - D f Since 
A m i n should scale as A (which we have also confirmed 
numerically); A m ; n ~ L~ z , we have 
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This generalizes the relation z = Df/{u) + 1) used as 
a criterion for a localized excitation by Lin et al. [2l| . 
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FIG. 4: (color online) (a) Distribution of the scaled gap ej = 
AiL a , with a — 2.84. The line corresponds to a gap exponent 
ui = 1. (b) Distribution of the upper bound A* scaled with 
dynamic exponent z — 3.6. The line shows the asymptotic 
small-gap behavior expected with uo = 1. The dashed curve 
is the Frechet distribution, Eq. 



With a = 2.84 and u> = 1 from the universal AiL a curve 
in the upper panel of Fig.0] we get z = 3.79, in excellent 
agreement with the value 3.77(4) obtained in Fig. [21 on 
the basis of (A*). The full probability distribution of 
A* is also very consistent with u) = 1, as shown in the 
lower panel of Fig.^ Here the best-fit dynamic exponent, 
z « 3.6, is slightly smaller than the value quoted above, 
but considering the limited number of cluster sizes and 
samples the agreement is quite satisfactory. We conclude 
that the scaling behavior is consistent with low-energy 
excitations due to finite local moment pairs. 

We can also calculate the full distribution of A. Using 



P M (e) = MP(e)[l- P< {e)] 



M-l 



P M (e) 



de 



3 -Mp<(e) 



(8) 



and Eq. © gives the Frechet distribution pfi|. 

P(u) = Au"exp(-A{uj + l)- 1 ^ 1 ), u = u AL z . (9) 

This form can only be fitted to the A* data in Fig. 01 
in the small-gap region. It is not completely clear to us 
why it fails for larger gaps, but there are signs of slower 
size-convergence on this side of the distribution. The 
deviations may also be related to the fact that we are 
studying the upper bound A*, not the exact gap A. 
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For clusters with ua 7^ rig, the final result of the 
RG moment decimation procedure would not be two mo- 
ments forming a singlet, but a single moment with spin 
S = \ua — ^b|/2. The global sublattice imbalance would 
necessarily lead to moments that on average grow as the 
decimation proceeds, unlike in the case of ha — tlb where 
it is possible that the effective moments reach a steady- 
state size as they are equally likely to grow or shrink at 
every RG step (a steady-state size is of course not guaran- 
teed on purely statistical grounds, but our scaling study 
above supports it). The lowest excitation for an tia =/= ng 
cluster must hence involve spins scattered through the 
whole cluster. Considering that (S) ~ ^/n, such an exci- 
tation cannot be local. We indeed find that the exponent 
relation J7J is violated for clusters with even n but no 
further restrictions on tia and Ub- We get a w 1.9 from 
scaling A^, which is considerably smaller than a 2.8 
for ha = 71b- However, the dynamic exponent obtained 
by scaling A* still remains z w 2Df. 

To conclude, we have presented a scaling study of the 
excitations of percolating 5 = 1/2 Heisenberg clusters, 
showing that the dynamic exponent z — 3.7(1); a factor 
2 larger than z — Df = 91/48 expected for quantum 
rotor states We propose that the low-energy excita- 
tions are due to interactions between effective moments 
formed by local imbalance in sublattice occupation. We 
have used concepts of strong-disorder RG and extremal 
value statistics to derive a scaling relation, Eq. J7J), be- 
tween z and two exponents governing the scaling of local 
gaps, under the condition that the low-energy excitations 
involve a finite number of spins. For clusters with an 
equal number of spins on the two sublattices (ua = tib), 



this relation is satisfied by numerical results for an up- 
per bound of the gap, A*, defined in Eq. @, and lo- 
cal gaps Ai similarly defined. However, for clusters with 
n-A 7^ n-s the global sublattice imbalance implies that the 
renormalized moments grow with system size. In accord 
with this non-locality we find that our exponent relation 
is violated. However, the dynamic exponent z w 2Df for 
both ua — tib and ua ^ Ub clusters. 

Our scenario implies that once two identical clusters 
are coupled in a bilayer, with a small inter-layer coupling 
J_l (smaller than the value at which the long-range order 
vanishes the low-energy excitations should change 

as effective moments in opposite layers comletely com- 
pensate each other. Preliminary calculations for clusters 
with J±/J ~ 0.05 indeed show a much smaller z, con- 
sistent with z — Df and the quantum rotor picture p|. 
This is also the case for the 3D classical Heisenberg model 
with columnar defects [2j| , which should describe the bi- 
layer but not the single layer for which there are uncom- 
pensated Berry phases (local sublattice imbalance) not 
captured by the mapping to the classical model. 

Another type of low-energy excitation possible in a per- 
colating cluster is a fracton 0, H3 • The fracton is lo- 
calized for finite energy f2, but delocalizes as — > 0. 
Hence, based on our scaling results for ua — tlb, it is 
different from our effective moment excitation. For a fi- 
nite system f2(£) oc L~ Df / df , where df is the fracton 
dimensionality. Approximate calculations of df pij have 
given Zf = Df/df smaller than the z reported here. 
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